Molecular surveillance of Plasmodium falciparum drug-resistance markers in Vietnam using multiplex amplicon sequencing (2000–2016)

Emergence and spread of Plasmodium falciparum resistance to artemisinin-based combination therapies (ACT) is a major challenge for Greater Mekong Subregion countries in their goal to eliminate malaria by 2030. Tools to efficiently monitor drug resistance beyond resource-demanding therapeutic efficacy studies are necessary. A custom multiplex amplicon sequencing assay based on Illumina technology was designed to target the marker of partial resistance to artemisinin (K13), five candidate modulators of artemisinin resistance, the marker of resistance to chloroquine (crt), and four neutral microsatellite loci. The assay was used to genotype 635 P. falciparum-positive blood samples collected across seven provinces of Vietnam and one of Cambodia between 2000 and 2016. Markers of resistance to artemisinin partner-drugs piperaquine (copy number of plasmepsin-2) and mefloquine (copy number of multidrug-resistance 1) were determined by qPCR. Parasite population structure was further assessed using a 101-SNP barcode. Validated mutations of artemisinin partial resistance in K13 were found in 48.1% of samples, first detection was in 2000, and by 2015 prevalence overcame > 50% in Central Highlands and Binh Phuoc province. K13-C580Y variant became predominant country-wide, quickly replacing an outbreak of K13-I543T in Central Highlands. Mutations in candidate artemisinin resistance modulator genes paralleled the trends of K13 mutants, whereas resistance to piperaquine and mefloquine remained low (≈ 10%) by 2015–2016. Genomic tools applied to malaria surveillance generate comprehensive information on dynamics of drug resistance and population structure and reflect drug efficacy profiles from in vivo studies.


Supplementary Methods
Design of the amplicon sequencing panel. SNP targets were selected for their reported association with antimalarial drug resistance in South-East Asia at time of experimental design (June 2017) and included: K13 as a marker of ART-R (full-length gene) 1 , exonuclease codon E415G as marker of PPQ resistance 2 , crt codons 72-76 as markers of chloroquine resistance 3 , mdr1 codon Y184F -the most frequent mdr1 SNP in SEA-as a modulator of resistance to multiple antimalarials 3,4 , as well as potential genetic modulators of MAL10:688956 and MAL13/RAD5-homolog-S1158A 5-8 . Genetic diversity was characterized using four microsatellite markers (poly-α, ARAII, TA81 and pk2) selected based on their high and exonuclease-E415G resulting in a total of 21 amplicons for 11 genes (see Supplementary Tables S1 and S2). The assay covered the whole propeller domain of K13 and 73% of the upstream sequence between codons 1-440 (see Supplementary Fig. S2). oligonucleotides program was run in a Veriti thermocycler (Applied Biosystems). After washing of unbound oligos and extension-ligation, products were indexed (i5 and i7) and PCR amplified (95°C for 3 min, followed by 33 cycles at 98°C for 20 sec, 67°C for 20 sec, 72°C for 40 sec, and final 72ºC for 1 min). A 5μl aliquot of purified PCR products was run in 2% agarose gel to confirm the correct amplicon size of 380-400bp. Each library was quantified using KAPA quantification kit for LightCycler 480 (Roche) and normalized to 5nM. All samples in a 96-well plate were pooled and stored at -20ºC. Prior to sequencing, pools were thawed, re-quantified using Qubit and adjusted to 4nM with Tris HCl according to the formula: (concentration in ng/μl / (660 g/mol x 380bp)) x 10 6 . Libraries were denatured with an equal volume of 0.2N NaOH to a loading concentration of 13pM.

Selective whole genome amplification (sWGA
Sequencing was conducted at Centre for Medical Genetics (University of Antwerp, Edegem, Belgium) in a MiSeq instrument with MiSeq Reagent Kit v2 (Illumina) and pairedend sequencing of 2x150bp reads. Sequence data analysis. Demultiplexing, alignment and variant calling were performed using TruSeq Amplicon Workflow in MiSeq Reporter software (Illumina). The genome sequence of P. falciparum 3D7 PlasmoDB build 29 was used as reference. Variants were identified using Genome Analysis Toolkit (GATK, Broad Institute) and a variant filter quality cut-off score of 30. Individual VCF files were combined with Picard tools (Broad Institute), and the allele frequency and depth information were summarized in Excel spreadsheets.
Genome Browse (Golden Helix) was used for visualization. Read depth cut-off was set at the maximum number of mapped reads for each amplicon found in results from 12 negative controls to control for spurious amplification. Mutations were only reported if found in at least two samples, or if the read count for the alternative allele was above the read cut-off for that amplicon. Loci with missing calls in more than 50% of all samples and controls were excluded. Haplotypes were built from calls with a ≥75% within-sample allele frequency to minimize risk of confounding by complex infections. Two amplicons (k13.i and MAL10) were sequenced in less than 50% of positive controls and were excluded from further analysis (see Supplementary Table S2). Median read depth of positive controls for the remaining amplicons in the assay was 7074x (medians range 34-27805x (Supplementary Table S2). Sensitivity for minor allele calling was estimated using three combinations of P. falciparum strains with known drug-resistance genotypes (3D7, wild-  The objective of this procedure is to clean up sWGA amplification products using magnetic beads. This step should be conducted in a separate area. -Bring AMPureXP beads, samples, and Resuspension buffer (Illumina) to room temperature (allow for at least 30 min after you take out of the fridge). TE buffer can be used alternatively to Resuspension buffer. -Prepare fresh 80% ethanol (400 μl/sample). o Allow reagents in bead cleaning procedures to reach room temperature for at least 30min. Vortex and mix well several times, including in-between pipetting steps.
o For all incubations leave the plate in the magnet for 5 min. Use an elastic rubber to hold the plate in magnet.
o All ethanol solutions should be freshly prepared.

DAY 3. Library Preparation (quantification, hybridization, extensionligation, amplification) -est. time full day
-Program the Hybridization program in the Veriti thermocycler before starting (please refer to Illumina's Reference Guide for programs in alternative thermocyclers) -Clean lab bench and all materials with a DNA decontamination solution (DNA Zap) and wipe out with water. -Identify the necessary reagents for hybridization step and leave them aside.
-Measure the ng/μl of purified products and other samples to be used in library preparation using Qubit HS DNA kit. o Targeted products should have a length of 380-400 bp.
-Store LNP plate at -20ºC until library quantification and sequencing (or proceed to next step). The library quantification using KAPA kit starts here. Procedures are based on KAPA Technical Data Sheet (reference list). Note that two plates will be necessary to quantify a full libraries plate. After loading the first qPCR plate, start with the second plate dilutions.
- -Take 5μl of each row using multichannel and transfer to a 8-row PCR strip. Then mix the contents of each well into one single tube (label "PAL-RoundX 5 nM"). -Store left-overs of normalized SGP plate at -20ºC. -Divide PAL tube in two aliquots and store at -20ºC, or proceed to the next step if sequencing is done in the next day. -Perform PAL concentration adjustment to 4nM as close as possible to the sequencing run.